#ifndef __fe_hh__
#define __fe_hh__

#include "femlib_config.hh"

#include "point.hh"
#include "fem_types.hh"

template<typename T>
struct FETraits
{
};

template <unsigned int DIM, fastfem::ElementType ELEMTYPE, unsigned int ORDER, fastfem::PolynomialFamily PFAMILY, fastfem::FunctionSpace FETYPE, typename T=double>
class FiniteElement
{
public:
  FiniteElement()  { }
  virtual ~FiniteElement() { }
  
  /** virtual function which evaluates the values of all
       shape functions at a given point ip */
  virtual T evaluate_basis (  unsigned int , const Point<SPATIAL_DIMENSION, T>& ip) = 0 ;
  
  /** virtual function which evaluates the values of all
       shape functions at all point ips */
  virtual void evaluate_basis_all ( const PointList<SPATIAL_DIMENSION, T>& ips, std::vector<T>& basis_vals ) = 0 ;
  
#ifdef FADBADPP_GRADIENTS  // see the note in fe.cc to enable implementation below
  virtual void evaluate_gradient ( unsigned int ib, const Point<SPATIAL_DIMENSION, T>& pt, Point<SPATIAL_DIMENSION, T>& grad) = 0 ;

  virtual void evaluate_gradient_all ( const PointList<SPATIAL_DIMENSION, T>& ips, std::vector<Point<SPATIAL_DIMENSION, T> >& grad_vals ) = 0;
#else
  virtual void evaluate_gradient ( unsigned int ib, const Point<SPATIAL_DIMENSION, T>& pt, Point<SPATIAL_DIMENSION, T>& grad) = 0 ;
  
  virtual void evaluate_gradient_all ( const PointList<SPATIAL_DIMENSION, T>& ips, std::vector<Point<SPATIAL_DIMENSION, T> >& grad_vals ) = 0 ;
#endif

#ifdef FADBADPP_HESSIANS  // see the note in fe.cc to enable implementation below
  virtual void evaluate_hessian ( unsigned int ib, const Point<SPATIAL_DIMENSION, T>& pt, PointList<SPATIAL_DIMENSION, T>& hess) = 0 ;

  virtual void evaluate_hessian_all ( const PointList<SPATIAL_DIMENSION, T>& ips, std::vector<PointList<SPATIAL_DIMENSION, T> >& hess_vals ) = 0;
#else
  virtual void evaluate_hessian ( unsigned int ib, const Point<SPATIAL_DIMENSION, T>& pt, PointList<SPATIAL_DIMENSION, T>& hess) = 0 ;

  virtual void evaluate_hessian_all ( const PointList<SPATIAL_DIMENSION, T>& ips, std::vector<PointList<SPATIAL_DIMENSION, T> >& hess_vals ) = 0 ;
#endif
  virtual int test() = 0 ;
  
//  virtual double* reference_to_physical_coords ( int n_points, double ref[] ) ;
//  virtual double* physical_to_reference_coords ( int n_points, double phy[] ) ;
  
  virtual void reference_element_coords ( std::vector<Point<SPATIAL_DIMENSION, T> >& nodes ) = 0 ;


};


#endif // #ifndef __fe_hh__
